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Reactions  involving  N  and  O  atoms  dominate  the  energetics  of  the  reactive  air  flow  around  space¬ 
craft  when  reentering  the  atmosphere  in  the  hypersonic  flight  regime.  For  this  reason,  the  thermal 
rate  coefficients  for  reactive  processes  involving  0(^P)  and  NO(^n)  are  relevant  over  a  wide  range 
of  temperatures.  For  this  purpose,  a  potential  energy  surface  (PES)  for  the  ground  state  of  the  NO2 
molecule  is  constructed  based  on  high-level  ab  initio  calculations.  These  ab  initio  energies  are  rep¬ 
resented  using  the  reproducible  kernel  Hilbert  space  method  and  Legendre  polynomials.  The  global 
PES  of  NO2  in  the  ground  state  is  constructed  by  smoothly  connecting  the  surfaces  of  the  grids  of 
various  channels  around  the  equilibrium  NO2  geometry  by  a  distance-dependent  weighting  function. 
The  rate  coefficients  were  calculated  using  Monte  Carlo  integration.  The  results  indicate  that  at  high 
temperatures  only  the  lowest  A-symmetry  PES  is  relevant.  At  the  highest  temperatures  investigated 
(20  000  K),  the  rate  coefficient  for  the  “0102-t-N”  channel  becomes  comparable  (to  within  a  factor 
of  around  three)  to  the  rate  coefficient  of  the  oxygen  exchange  reaction.  A  state  resolved  analysis 
shows  that  the  smaller  the  vibrational  quantum  number  of  NO  in  the  reactants,  the  higher  the  rela¬ 
tive  translational  energy  required  to  open  it  and  conversely  with  higher  vibrational  quantum  number, 
less  translational  energy  is  required.  This  is  in  accordance  with  Polanyi’s  rules.  However,  the  oxygen 
exchange  channel  (N02-I-01)  is  accessible  at  any  collision  energy.  Einally,  this  work  introduces  an 
efficient  computational  protocol  for  the  investigation  of  three-atom  collisions  in  general.  ©  2014  AIP 
Publishing  LLC.  [http://dx.doi.Org/10.1063/l.4897263] 

I.  INTRODUCTION 

Nitric  dioxide  (NO2)  plays  a  major  role  in  atmo¬ 
spheric  chemistry,  as  a  smog  constituent  and  in  combustion 
processes.'"^  In  the  atmosphere,  nitrous  acid  (HONO)  can  be 
generated  from  NO2  through  reaction  with  water  which  is  im¬ 
portant  because  HONO  is  a  major  source  of  OH  leading  to 
the  formation  of  ozone  and  other  air  pollutants."^’ ^  Depending 
on  the  chemistry  involved,  HONO  formation  from  NO2  can 
also  involve  reactants  such  as  NO,  water  on  surfaces^  or  soot 
surfaces.®  While  most  of  these  processes  occur  close  to  ther¬ 
mal  equilibrium  and  at  moderate  temperatures,  reactions  in¬ 
volving  NO2  have  also  attracted  interest  under  more  extreme 
conditions.  This  is  the  case  for  the  hypersonic  flight  regime 
of  spacecraft  reentering  the  atmosphere.  The  chemistry  near 
the  surface  of  such  vehicles  typically  involves  highly  non¬ 
equilibrium  conditions  with  vibrational  and  rotational  temper¬ 
atures  reaching  several  thousand  Kelvin.^  The  gas-phase  and 
surface  reactions  and  energy  transfer  at  these  temperatures  are 
essentially  uncharacterized  and  the  experimental  methodolo¬ 
gies  capable  of  probing  them  are  not  well  established.  Un- 
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der  such  circumstances,  validated  computational  investiga¬ 
tions  become  a  valuable  complementary  tool. 

In  the  present  work,  the  formation  of  NO2  from  NO(^n) 
and  0(^P)  and  its  subsequent  dissociation  into  oxygen- 
exchange  and  O2  formation  channels  is  investigated.  Char¬ 
acterizing  the  kinetics  and  dynamics  of  a  chemical  reaction 
from  computation  usually  requires  a  potential  energy  surface 
(PES)  for  the  interactions  involved.  Within  the  framework  of 
Molecular  Dynamics  (MD)  simulations,  this  can  either  be  an 
ab  mfiio/semiempirical  MD  (AIMD)  approach,  a  model  PES 
or  a  parameterized,  multidimensional  PES  fitted  to  electronic 
structure  calculations,  experimental  observables,  or  a  com¬ 
bination  of  the  two.®  Using  AIMD  is  usually  only  possible 
for  small  systems  and  sufficiently  short  simulation  times  due 
to  the  considerable  computational  cost.'^ Eor  parametrized 
PESs,  the  selection  of  the  grid  points  on  which  the  ab  initio 
calculations  are  carried  out  is  always  a  compromise  between 
the  computational  cost  and  the  intended  purpose  of  covering 
the  most  significant  regions  of  the  configuration  space  (which 
is  a  priori  unknown).  A  more  efficient  approach  is  to  use  sep¬ 
arate  grids  for  the  different  regions  of  interest  and  join  them 
smoothly.  In  the  asymptotic  regions,  that  is,  reactants  or  prod¬ 
ucts,  the  dynamics  is  carried  out  on  a  single  and  simpler  PES 
leading  to  a  faster  simulation.  In  the  intermediate  regions,  the 
global  surface  is  represented  as  a  weighted  sum  of  the  contri¬ 
butions  of  all  individual  PESs.  This  approach  has  been  suc¬ 
cessfully  applied  to  connect  different  force  fields.'^  We  will 
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pursue  a  similar  strategy  in  the  present  work  to  build  a  global 
PESforN02. 

A  global  PES  for  the  l^Aj  and  1^62  electronic  states 
of  NO2  based  on  high-level  ab  initio  calculations  has  been 
computed  some  years  ago.'^  Very  good  agreement  with  the 
experimental  results  for  the  dissociation  energies'^  and  the 
equilibrium  geometries  was  obtained  from  this  PES.^^  The  re¬ 
sulting  frequencies  calculated  using  filter  diagonalization  rea¬ 
sonably  matched  the  experiments  and  values  simulated  from 
other  published  PESs.  Both  PESs  in  Ref.  13  were  represented 
with  three-dimensional  cubic  splines.  In  subsequent  work,  the 
long-range  part  of  the  l^Aj  PES  and  the  lowest  l^A"  PES 
were  corrected  by  considering  the  NO  potential  energy  curve 
as  a  separate  term  in  the  PES  at  large  NO-0  separations.  With 
this  modification,  the  modeling  of  the  collision  dynamics  to 
study  the  exchange  and  recombination  reactions  in  O-I-NO 
collisions  was  performed.'®  For  temperatures  between  100 
and  2500  K,  the  recombination  (kJo)  and  exchange  rate  co¬ 
efficients  (kS)  varied  slightly  within  a  factor  of  2.*^ 

The  aim  of  the  present  work  is  the  characterization  of  the 
reaction  between  0(^P)  and  NO(^n)  at  higher  temperatures 
relevant  to  the  hypersonic  flight  regime  of  reentering  space¬ 
crafts.  At  a  more  fundamental  level,  we  are  particularly  in¬ 
terested  in  establishing  an  efficient  computational  framework 
capable  of  probing  small  molecular  systems  for  a  wide  range 
of  parameters.  The  general  approach  is:  (1)  the  construction 
of  PESs  for  the  ground  state  of  the  NO2  molecule  based  on 
high-level  ab  initio  calculations  and  their  representation  with 
a  reproducing  kernel  Hilbert  space  (RKHS)  method  combined 
with  Legendre  polynomials;  (2)  quasi-classical  trajectory  cal¬ 
culations  to  study  the  adiabatic  reaction  dynamics,  and  (3) 
calculation  of  the  rate  coefficients  for  the  different  exit  chan¬ 
nels  using  a  Monte  Carlo  method.  The  results  of  the  sim¬ 
ulations  will  be  compared  with  previous  computations'®""' 
and  experiments'^  in  order  to  validate  the  computational  ap¬ 
proach.  In  the  context  of  the  hypersonic  flight  regime,  the  O2 
production  channel  (NO-fO— >-02-|-N)  is  of  particular  interest 
since  the  thermal  dissociation  threshold  for  O2  is  lower  than 
that  for  N2  in  air.  A  central  question  is  how  translational  ex¬ 
citation  of  the  reactants  and  vibrational  excitations  of  the  NO 
affect  the  rate  coefficient  of  O2  formation. 

A  major  focus  here  is  to  test  and  implement  a  compu¬ 
tational  protocol  which  should  be  suitable  for  routine  and 
largely  automated  procedures  in  determining  important  fun¬ 
damental  inputs  for  reaction  networks  relevant  to  atmospheric 
and  astrophysical  modeling.  The  characterization  of  more 
complex  chemical  systems  involves  the  modeling  of  a  collec¬ 
tion  of  reactions  which  can  be  arranged  in  a  reaction  mech¬ 
anism  where  multiple  elementary  reactions  take  place  at  a 
wide  range  of  changing  conditions  such  as  temperature  and 
pressure. In  most  cases,  chemical  processes  in  these  net¬ 
works  occur  out  of  equilibrium.  This  complicated  panorama 
makes  difficult  the  possibility  of  carrying  out  experiments 
nowadays.  Even  in  implementing  large  scale  simulations,  a 
significant  amount  of  data  as  input  is  required,  for  instance, 
the  rate  coefficients  of  the  elementary  reactions.  This  infor¬ 
mation  may  come  from  different  sources  either  theoretical  or 
experimental  and  it  happens  frequently  that  important  data 
are  not  available.  Given  this  scenario,  the  search  for  efficient 


computational  frameworks  for  the  calculation  of  the  rate  co¬ 
efficients  is  an  issue  of  primary  interest. 

The  outline  of  the  present  work  is  as  follows.  In  Sec.  II, 
the  methodology  for  constructing  the  PESs,  the  trajectory  cal¬ 
culations,  and  the  computation  of  the  thermal  quantities  is 
discussed.  Section  III  presents  and  discusses  the  results  and 
in  Sec.  IV,  conclusions  about  the  validity  of  the  approach  are 
drawn. 

II.  METHODOLOGY 
A.  Ab  initio  calculations 

The  present  calculations  use  Dunning’s  standard 
correlation-consistent  polarized  quadruple  (cc-pVQZ)  basis 
set,  which  is  sufficiently  large  to  properly  describe  the 
N-O  and  0-0  bonds. ^®’^®  For  a  meaningful  description 
of  the  geometries  and  energetics  of  the  NO2  molecule,  its 
dissociation  energetics  and  the  conical  intersections  between 
the  ground  and  the  low-lying  excited  states,  a  multiconfig- 
urational  method  is  required.  For  this,  the  complete  active 
space  self-consistent  field  (CASSCF)^^  method  is  employed 
which  takes  into  account  static  correlation.  For  recovering 
the  dynamic  correlation,  e.g.,  during  bond  breaking  and 
formation,  multireference  configuration  interaction  (MRCI) 
is  used.  The  Davidson  correction^^’^®  is  included  to  estimate 
the  effects  of  higher  excitations  (MRCI-fQ).  All  calculations 
are  performed  with  the  MOLPRO  program  ver.  2012.1.^® 

Figure  1  gives  a  schematic  representation  of  the  general 
relevant  stages  of  the  reaction  of  interest.  NO2  formation  has 
a  barrier  of  1.3  eV  from  the  N-I-0102  side,  whereas  from  the 
NOl-l-02  side  it  is  a  barrierless  process.  Initially,  the  0(®P) 
atom  approaches  the  NO(^n)  molecule.  Depending  on  the 
impact  parameter  b,  the  quantum  state  of  the  NO  (vibrational 
and  rotational)  and  the  collision  energy,  the  NO2  molecule  is 
formed.  In  the  absence  of  additional  collisional  partners  (zero 
pressure),  the  NO2  will  decay  as  it  was  formed  with  an  energy 
above  its  dissociation  limit.  The  following  possibilities  exist 
for  the  products:  NO2  can  (1)  dissociate  into  the  same  incom¬ 
ing  constituents  (NOl-l-02  ^  NOl-l-02),  (2)  dissociate  into 
the  same  constituents  but  with  the  oxygen  atoms  exchanged 
(N01-f02  ^  N02-f01),  (3)  form  dioxygen  (O2+N),  or  (4) 
atomize  into  the  individual  atoms  (N-|-01-|-02). 

In  the  presence  of  additional  collisional  partners,  the  NO2 
can  loose  its  excess  energy  and  stabilize.  In  the  high-pressure 


O2+N  ^  NO2  ^  NO  +  O 


FIG.  1.  Relevant  stages  in  the  NOl+02  reaction.  The  energies  were  derived 
from  the  ab  initio  calculations  of  this  work. 
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limit,  this  occurs  for  all  trajectories  which  enter  the  potential 
well  of  NO2  and  thereby  it  is  straightforward  to  determine  the 
high-pressure  limit  rate  coefficient  for  NO2  formation. 

In  order  to  describe  the  different  channels,  three  grids  in 
Jacobi  coordinates  are  employed  to  systematically  cover  the 
relevant  regions  of  the  PES  for  NO2  (Figure  1).  The  hrst  grid 
(for  the  NO  1+02  asymptote)  includes  the  coordinates  (the 
N-01  distance),  (the  distance  from  the  center  of  mass  of 
NOl  to  the  oxygen  atom  02),  and  (the  angle  formed  by 
01,  the  center  of  mass  of  NOl  and  02  (Figure  1)).  With  this 
grid,  the  entrance  (N01+02)  and  N02-formation  channels 
can  be  described.  For  the  oxygen  exchange  channel  (N02- 
01),  a  permutation  of  the  oxygens  in  the  first  grid  is  per¬ 
formed.  As  a  result,  the  second  grid  is  dehned  by  ^2,  /?2,  and 
cTj-  For  the  third  grid  (0102+N),  r3  is  the  0-0  distance, 
is  the  distance  from  the  center  of  mass  of  the  O2  molecule  to 
the  nitrogen  atom,  and  Q'3  is  the  angle  between  02,  the  center 
of  mass  of  the  O2  molecule,  and  the  nitrogen  atom  (Figure  1). 

For  the  and  r3  coordinates,  12  points  (n^  =12)  from 
1.7  to  2.55  flQ  were  selected.  These  points  include  the  equi¬ 
librium  positions  and  the  turning  points  of  the  hrst  four  vi¬ 
brational  levels  of  the  NO  and  O2  molecules,  respectively. 
Three  additional  points  were  used  for  a  better  description  of 
the  ground-state  PES  wells  of  the  NO  and  O2  molecules.  The 
turning  points  for  the  NO  molecule  were  those  of  a  spec¬ 
troscopically  accurate  Rydberg-Klein-Rees  potential.^^’^^  For 
the  O2  molecule,  the  turning  points  were  obtained  from  the 
solution  of  the  Schrodinger  equation  for  a  Morse  poten¬ 
tial  htted  to  density  functional  theory  data.^^  As  for  the  Ry 
and  R3  coordinates,  29  points  (n^j  =  29)  between  1.85  and 
18.9  Aq  were  selected.  The  and  Q'3  angles  (in  degrees) 
are  from  an  11 -point  (n„  =  11)  Gauss-Legendre  quadrature 
(a  =  11.98,  27.49,  43.10,  58.73,  74.36,  90.0,  105.64,  121.28, 
136.9,  152.5,  168.02).  In  total,  3828  points  were  calculated 
for  each  grid. 

B.  Representation  of  the  PESs 

In  order  to  obtain  a  continuous  PES  for  the  NO2  molecule 
from  the  reference  ab  initio  calculations,  a  RKHS  ansatz^'^ 
(see  below)  for  the  distances  (r  and  R)  is  combined  with  Leg¬ 
endre  polynomials  for  the  angles  (a)  (see  Figure  1).  First,  for 
each  combination  (R,a),  a  set  of  functions  ^{r)  is  obtained 
by  using  a  RKHS.  In  a  second  step,  since  the  ab  initio  calcula¬ 
tions  were  carried  out  at  Gauss-Legendre  points,  the  2d-PES 
W{R,a)  can  be  expanded  in  terms  of  Legendre  polynomials 
(cos  O')  as  follows: 

10 

W(R,a)^'^f)^(R)P^(cosa),  (1) 

;.=o 

where  fx(R)  are  radial  coefficients  (“radial  strength 
functions”)^^  and  X  is  their  order.  A  second  RKHS  in¬ 
terpolation  is  then  used  along  the  P-coordinate. 

The  energy  at  an  off-grid  point  (/,  R',  and  a')  is  evaluated 
as  follows.  From  the  direct  evaluation  of  the  functions  Vj^  „(r) 
at  F,  a  set  of  energies,  W(R,  a',  r')  for  each  combination  of 
(R,  a)  on  the  grid  is  obtained  which  yields  a  2-dimensional 
PES.  By  expanding  W(R,  a',r')  in  Legendre  polynomials  up 
to  order  =10,  the  following  n^  x  n^  linear  system  is 


obtained: 

1 

^  =  'F(P^-,  0(^;r'),  (2) 

;i=o 

where  q;^.,  k—  1,  ■  ■  ■ ,  {n^  —  !)>  are  the  Gauss-Legendre  points. 
Solving  this  linear  system  for  f^{Rp  on  all  grid  points  R- 
(see  Eq.  (2))  yields  new  radial  strength  functions  which  are 
represented  as  a  RKHS.  These  kernels  can  be  evaluated  at 
position  R' .  Finally,  evaluating  Eq.  (2)  at  =  a'  gives  the 
energy  W{R' ,  a'-,  r')  at  an  off-grid  point.  It  should  be  empha¬ 
sized  that  truncating  the  Legendre  expansion  at  hnite 
potentially  leads  to  inaccuracies  in  regions  where  variations 
along  a  are  appreciable.  However,  it  is  found  that  radial  cuts 
along  R  for  neighboring  angles  show  similar  behaviours. 
An  explicit  comparison  is  provided  in  the  inset  of  Figure  4. 
Based  on  this,  retaining  1 1  terms  in  the  angular  degree  of  free¬ 
dom  was  deemed  sufficient. 

Details  about  the  RKHS  procedure  can  be  found  in  the 
literature.^"^  Here,  only  the  essentials  are  summarized  for  com¬ 
pleteness.  Within  the  RKHS  framework,  a  multidimensional 
function  F(R)  (here  the  PES)  which  depends  on  the  molec¬ 
ular  coordinates  R  is  assumed  to  belong  to  a  RKHS,  that  is, 
V (R)  is  a  bounded  linear  function.  Thus,  the  function  V (R) 
can  be  represented  as 

M 

F(R)  =  ^«,e(R„R),  (3) 

k=l 

where  R  =  (xj,  . . . ,  x'p)  indicates  the  F  internal  coordinates 
for  a  given  nuclear  conhguration  of  the  system,  at  which  the 
interpolated  energy  is  sought,  R^.  =  (x^j ,  . . . ,  x^^^)  are  the  co¬ 
ordinates  of  the  M  ab  initio  energies,  and  2(Rjfc>  R)  is  the  re¬ 
producing  kernel  which  consists  of  a  set  of  linearly  indepen¬ 
dent  functions  also  belonging  to  the  RKHS.  The  in  Eq.  (3) 
satisfy  the  relation 

M 

=  y(R,),  k=l,2,...,M,  (4) 

l=\ 

where  =  g(R^,  R^)  and  V (R^)  are  the  ab  initio  energies  at 
the  grid  points  (Figure  1).  The  a,  are  determined  once  from  a 
Cholesky  decomposition.^®  The  reproducing  kernel,  QiR^,  R) 
for  a  multidimensional  PES  can  be  represented  as  a  product 
over  1 -dimensional  kernels  x[) 

e(R„R)=  nq,{x,,,x[).  (5) 

i=\ 

For  1 -dimensional  distance-like  variables,  the  following  ex¬ 
plicit  expression  for  the  kernel  was  employed: 

g"’'”(x,  x')  —  B(m  +  1,  n) 

X  2^3  M  +  1,  nr  +  l;n  +  m  +  1;  — (6) 

where  x^  and  x^  are  the  larger  and  the  smaller  of  x  and 
xf,  respectively,  B{m  +  1,  n)  is  the  beta  function,  and 
2+1  (—n  +  1,  OT  +  1;  n  +  OT  +  1;  ^)  is  the  Gauss  hypergeo¬ 
metric  function. The  superscripts  n  and  m  are  related  to  the 
smoothness  of  the  function  to  be  interpolated  and  the  recip¬ 
rocal  power  of  the  distance  weighting  factor,  respectively.  In 
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TABLE  I.  Spectroscopic  constants  for  the  ground  state  of  NO  and  O2 
molecules. 


NO  molecule 

Our  calculations 

Experiments  (Ref.  33) 

(cm“^) 

1922.01 

1904.1 

(cm“*) 

14.48 

14.09 

O2  molecule 
(cm~^) 

1607.5 

1580.2 

(cm“^) 

14.1 

12.0 

the  present  work,  n  —  2  (smoothness  criterion  applied  up  to 
the  first  derivatives  of  the  potential)  and  m  —  6  which  cor¬ 
responds  to  a  1/x®  long-range  interaction,  although  different 
dependencies  could  also  be  considered. 

The  RKHS  procedure  is  a  generic  and  parameter-free 
method.  It  can  be  used  for  constructing  n-dimensional  PESs 
from  either  regular  or  irregular  grids  which  is  an  additional 
advantage.  The  asymptotic  behavior  of  the  interactions  can 
be  explicitly  taken  into  account  (see  above).  If  molecular 
symmetry  is  present,  this  can  also  be  employed  (e.g.,  per- 
mutational  symmetry)  and  explicit  analytical  derivatives  of 
the  kernel  are  available.  Because  the  long-range  behavior 
can  be  controlled  and  analytical  derivatives  are  available,  a 
RKHS-based  PES  offers  potential  advantages  over  cubic  (or 
higher-order)  spline  interpolations.  A  possible  drawback  of 
the  RKHS  method  is  the  cost  of  the  calculation  of  the  in 
Eq.  (4)  which,  however,  needs  to  be  done  only  once  and  the 
evaluation  of  the  reproducing  kernel  (2(Ri:,  R)  at  every  new 
calculation  of  the  potential  energy  V (R).  This  may  not  be  a 
problem  if  the  pre-summations  of  the  method  are  handled  in 
an  efficient  manner. 

To  properly  describe  the  asymptotic  behavior  of  the  PES 
for  large  separations  between  the  reactants,  Morse  potentials 
(V  =  D^[\  —  exp(— /8(r  —  were  fitted  to  the  ab  initio 

energies  for  the  NO  and  O2  monomers.  The  resulting  param¬ 
eters  were  —  7.91  eV,  p  —  1.34  Oq  \  and  —  2.17  Oq 
for  NO  and  =  5.67  eV,  P  =  1.37  and  =  2.29 
for  O2.  The  reported  spectroscopic  data  are  —  6.623  eV 
and  =  2.18  Oq  for  NO,  and  5.21  eV  and  2.28  Gq  for  O2, 
respectively.^^  The  disagreement  in  the  dissociation  energies 
of  NO  (experimental  and  calculated)  does  not  affect  the  re¬ 
sults  since  the  full  atomization  process  is  not  considered  here 
as  explained  further  below.  Eor  further  characterization  of  the 
NO  and  O2  molecules,  the  spectroscopic  constants  calculated 
from  the  Morse  potential  parameters  are  compared  to  experi¬ 
ments  in  Table  I.^°  The  results  indicate  a  better  agreement  for 
the  anharmonic  term  (cu^x^)  for  the  NO  molecule. 

For  the  global  PES  for  the  system,  the  local  PESs  deter¬ 
mined  for  the  three  grids  corresponding  to  the  relevant  chan¬ 
nels  of  the  reaction  (Figure  1)  need  to  be  joined.  For  this,  the 
global  PES  is  written  as  the  weighted  sum  of  the  three  asymp¬ 
totic  PESs 

3 

P(x)  =  ^«;/x)y,.(x),  (7) 

J=i 

where  x  =  (r[2,  r-^^,  r2^)  are  the  three  interatomic  distances 
of  the  NO2  molecule  (Figure  1),  and  ^[(x),  V'2(*)’  Vj{x) 


are  the  potential  energies  for  the  three  regions  for  which  the 
distances  r[2  (N-01),  (N-02),  and  r23  (01-02)  are  small, 

and  the  remaining  distances  vary  from  small  to  large,  respec¬ 
tively.  The  weights  uij{x)  are  coordinate-dependent  and  result 
from  the  normalization  of  the  following  weight  functions: 

Wjoixj)  =  exp(-x^./AR),  (8) 


E  w,oixi) 

/=i 

where  Xj  are  either  r[2  or  rj3  or  ^23  and  AR  =  0.05  Gq.  A  sim¬ 
ilar  weighting  function  was  successfully  used  in  Multisurface 
Adiabatic  Reactive  Molecular  Dynamics  (MS-ARMD)  which 
uses  the  energy  as  the  control  parameter  in  connecting  multi¬ 
ple  force  fields.'^ 

C.  Molecular  dynamics  simulations 
1.  Generation  of  the  initial  conditions 

Armed  with  a  globally  valid  PES,  the  dynamical  evo¬ 
lution  of  the  system  is  followed  by  propagating  Hamilton’s 
equations  of  motions  subject  to  initial  conditions.  The  equa¬ 
tions  of  motion  are  integrated  numerically  using  the  Velocity- 
Verlet  algorithm.^^  Since  the  total  linear  and  angular  momenta 
are  conserved  during  the  dynamics,  it  is  convenient  to  se¬ 
lect  the  center  of  mass  of  the  NO2  molecule  as  the  origin  of 
the  coordinate  system.  In  order  to  capture  the  fastest  nuclear 
movements  in  the  system,  a  time  step  of  2  a.u.  (a.u.  means 
atomic  unit  of  time)  (4.84  x  10“^  fs)  was  used  for  the  propa¬ 
gation.  Trajectories  were  started  from  an  initial  separation  of 
(min(15  Og,  impact  parameter  -|-  5  Gq))  and  followed  until  the 
fragments  had  separated  to  the  same  distances.  If  this  criterion 
was  not  met,  integration  continued  for  a  maximum  time  of  the 
estimated  interactionless  time  of  flight  (see  below)  through 
the  “reaction-zone”  plus  an  additional  100  ps,  during  which 
most  of  the  NO2  molecules  decayed  except  for  very  few  low- 
energy  (low  and  low  rovibrational  state  of  initial  NO) 
trajectories. 

Suitable  initial  conditions  for  NO  were  generated  from  a 
WKB'^^-quantized  periodic  orbit  of  the  corresponding  rotat¬ 
ing  Morse  oscillator  for  given  vibrational  v  and  rotational  j 
quantum  numbers."^^  The  symmetry  axis  of  the  NO  molecule 
and  the  axis  of  its  rotation,  which  are  orthogonal  to  each  other, 
was  also  randomly  orientated  as  was  the  NO  angular  momen¬ 
tum.  Sampling  of  the  initial  rovibrational  quantum  states  of 
NO,  the  impact  parameter  (b)  and  the  collision  energy  {EJ 
are  discussed  in  Sec.  II  D. 

For  characterizing  the  reaction  at  a  given  temperature,  4 
different  cases  from  I  to  IV  are  considered: 

(I)  corresponds  to  N02-formation  in  the  high-pressure 
limit.  The  NO2  molecule  stabilizes  through  collisions  with 
the  environment.  We  will  indicate  this  by  “00”  as  superscript 
(k“(T)).  For  the  formation  of  NO2,  the  same  criterion  as  in 
Ref.  16  is  assumed,  that  is,  the  oxygen  atom  02  is  captured 
once  it  approaches  to  within  less  than  3.78  Gq  of  the  center  of 
mass  of  NOl. 

(II)  is  for  NO2  formation  and  survival  which  carries 
the  label  “NO2”  as  superscript  (k^°2(7’)).  For  this  case,  the 
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“oxygen  capture  criterion”  (see  item  (I))  is  combined  with  a 
lifetime  criterion,  defined  as  >  0.  Here,  is 

the  total  time  of  the  trajectory,  which  is  either  (a)  —  100  ps 

for  the  case  that  NOj  does  not  dissociate  or  (b)  the  time  un¬ 
til  any  of  the  interatomic  separations  is  larger  than  the  ini¬ 
tial  separation  of  NO  1-02.  The  interactionless  time  of  flight 

hot — is  determined  from  the  initial  condi¬ 
tions. 

(III)  is  the  oxygen  atom  exchange  reaction,  “ex”  with  as¬ 
sociated  rate  coefficient  ^’‘(Tj.  The  criterion  for  this  case  is 
that  the  N-02  distance  is  small  and  the  N-01  or  01-02  dis¬ 
tances  are  large. 

(IV)  refers  to  oxygen  molecule  formation  (0102-l-N) 
and  the  corresponding  rate  coefficient  is  ^‘^2+^(7’).  This  case 
is  encountered  if  the  01-02  distance  is  small  and  the  N-Ol 
or  N-02  distances  are  large. 

An  additional  possible  process  is  full  atomization  into 
N-PO-fO.  However,  this  channel  is  not  considered  because 
even  at  T  =  20000  K,  when  the  NOj  complex  potentially 
contains  sufficient  energy  for  atomization,  the  probability  to 
decay  through  this  channel  is  small.  This  is  based  on  consid¬ 
ering  the  probability  distributions  for  the  rovibrational  energy 
of  the  NO  and  O2  molecules  in  the  product  channels  and  the 
following  reasoning. 

First,  direct  three-body  breakup  is  extremely  unlikely 
which  only  leaves  dissociation  into  an  atom  plus  diatom  sys¬ 
tem.  After  breakup  a  considerable  amount  of  energy  is  in  the 
relative  translation  of  the  atom  -p  diatom  products.  Dissocia¬ 
tion  of  the  diatomic  occurs  whenever  the  vibrational  energy 
(on  the  effective  potential)  exceeds  the  dissociation  energy 
level  of  the  diatomic  molecule.  When  the  rotational  quantum 
number  j  of  the  product  diatomics  is  smaller  than  j  =  50  (see 
inset  in  Figure  2),  the  effect  of  rotation  on  decreasing  the  dis- 


FIG.  2.  Probability  distributions  of  the  rovibrational  energy  for  NO  and  Oj 
molecules  in  the  products  at  20  000  K  are  shown  separately.  The  blue  and 
red  vertical  dashed  lines  are  the  calculated  dissociation  energies  for  O2  and 
NO,  respectively.  The  corresponding  experimental  values  are  6.62  eV  and 
5.21  eV  for  NO  and  O2,  respectively.  The  inset  in  the  hgure  shows  the  effec¬ 
tive  potential  of  the  O2  molecule  where  /r  is 

the  reduced  mass  between  the  oxygen  atoms  and  r  the  interatomic  distance) 
for  j  =  0,  50,  100,  150,  and  200. 


sociation  energy  via  centrifugal  forces  will  be  insignificant. 
Figure  2  reports  the  rovibrational  energy  distributions  for  NO 
and  02-  They  indicate  that  the  probability  for  NO-dissociation 
is  very  small  whereas  dissociation  of  O2  is  more  likely  (about 
20%).  However,  this  occurs  still  with  very  low  probability 
since  the  centrifugal  barrier  of  the  rotating  diatomics  does 
not  allow  the  molecule  to  dissociate  in  a  classical  simulation, 
even  if  its  energy  is  larger  than  the  asymptotic  energy  level 
for  decomposition  but  still  less  than  the  top  of  the  centrifugal 
barrier.  This  explains  why  fewer  dissociation  events  are  found 
than  expected  based  on  energetics  alone  and  is  confirmed  in 
the  simulation.  Only  63  out  of  the  10  000  trajectories  end  up 
in  N-pO-pO  and  all  of  them  take  place  via  O2  dissociation. 
The  implication  is  that  dissociation  into  either  the  O2  or  NO 
channels  is  fast  and  the  final  fragmentation  to  the  fully  atomic 
channel  is  slow.  Therefore,  the  atomic  channel  will  not  affect 
the  rate  coefficients  of  the  diatomic  channels. 


D.  Monte  Carlo  calculation  of  the  thermal  rate 
coefficients 

The  thermal  rate  coefficient'*^^*’  can  be  determined  from 

I -  00 

kiT„  7;)  = 

g{T^)y  nix  J 

0 

(10) 

where  =  kgTj  and  kg  is  the  Boltzmann  constant,  7j,  Tj.^, 
Tg  are  the  translational  temperature  of  NO  and  O,  the  rovi¬ 
brational  temperature  of  NO  and  the  electronic  tempera¬ 
ture  of  NO,  respectively,  g(T^)  is  the  electronic  degeneracy 
factor,  IX  is  the  reduced  mass  of  NOl  and  O2,  respec¬ 
tively,  and  o-(£(,;  Tj.^)  is  the  integral  cross  section  as  a  function 
of  the  collision  energy,  E^  and  In  our  studies,  it  was  as¬ 
sumed  that  Tj  =  Tj.^,  =  Tg.  In  other  words,  the  various  degrees 
of  freedom  in  the  reactants  are  in  thermal  equilibrium.  If  no 
subscript  is  shown  for  T,  then  it  corresponds  to  the  common 
temperature.  In  a  similar  manner,  the  reaction  cross  section 
cr(7?g;  Tj^)  for  a  given  collision  energy  E^  is 


X!  X!  vj) 

v=o  j=0 _ 

"max  4ax(") 

^  X!  +  l)e"^"'^”2 

v=0  j=0 


(11) 


where  <y^j{E^)  is  the  (v,  7)-state  dependent  cross  section  at 
collision  energy  E^.  The  energy  E^j  of  a  rovibrational  state 
(v,  J)  is  calculated  according  to  a  Morse  oscillator  model 
for  the  NO  molecule."*^  The  cross  section  as  an  integral  of 
the  opacity  function  P^j(b',  E^)  for  given  E^  and  rovibrational 
state  (v,  j)  is 


=  j  P,j{b-E^)2nbdb.  (12) 

0 

Assuming  a  thermal  distribution  of  the  reactants  0(^P) 
and  NO(^n)  the  expression  for  the  electronic  degeneracy 
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factor,  is 

g(T^)  =  [5  +  3  exp(-227.8  K/  T^)  +  exp(-326.6  K/  TJ] 

X  [1-f  exp(-177.1K/7;)].  (13) 


The  integral  in  Eq.  (10)  can  be  calculated  by  using  an 
Importance  Sampling  Monte  Carlo  scheme.^^  For  this,  the  vi¬ 
brational  and  rotational  quantum  numbers  v  and  j,  and  the  col¬ 
lision  energy  (E^)  are  sampled  from  the  following  probability 
distributions: 


(2 j  +  l)e 


Pvj(T,y)  ^ 

max  •'max'  ' 

E  E(2/+1K 

i>'=0  j'=0 


(14) 


4  E 

'^rv  v’ j’ 


and 


p(£,)d£,  =  (15) 


respectively. 

For  the  determination  of  the  upper  limit  (fe^ax)  of  the 
relevant  sampling  interval  for  the  impact  parameter,  b,  pre¬ 
liminary  tests  were  carried  out  at  low  temperatures  (100  and 
200  K).  As  NO2  was  not  formed  for  b  >  26  Uq,  the  impact 
parameter  was  uniformly  sampled  from  0  up  to  b^^^^  —  26  a^. 

After  evaluating  the  integral  in  Eq.  (10)  over  v,j,  and  E^, 
the  resulting  expression  for  k(T)  is 


k(T)^ 


I  8 

jtixp  g{T)N„ 


i=l 


(16) 


where  and  are  the  number  of  reactive  and  the  to¬ 
tal  number  of  trajectories,  respectively,  and  bj  is  the  impact 
parameter  of  reactive  trajectory  i.  The  convergence  of  the  in¬ 
tegral  in  Eq.  (10)  was  monitored  by  the  decrease  of  the  Monte 
Carlo  error. 


FIG.  3.  Contour  plots  for  the  NOl+02  (panel  (a))  and  0201+N  (panel  (b)) 
grids  in  R^,  and  Jacobi  coordinates,  respectively.  The  N-01  (rj) 

and  02-01  (r3)  distances  were  set  to  2.25  a^.  In  both  cases,  the  contour  levels 
are  separated  by  0.5  eV. 


Ml.  RESULTS  AND  DISCUSSION 
A.  The  ground  state  1^A'  PES 

In  the  present  work,  the  reaction  dynamics  of  NO2  on  the 
ground-state  (l^A')  PES  is  considered.  Even  when  the  excited 
states  can  approach  the  ground  state  PES  for  certain  config¬ 
urations  (see  Figure  1  in  the  supplementary  material"*®),  their 
influence  can  be  considered  to  be  minor  as  it  requires  6  eV 
to  access  them.  According  to  the  Maxwell-Boltzmann  dis¬ 
tribution  for  20  000  K,  translational  energies  beyond  6  eV 
are  accessed  with  little  probability  (0.06).  This  is  further  de¬ 
creased  if  one  takes  into  account  that  some  energy  goes  into 
rotational  and  vibrational  states  of  NO  which  are  distributed 
according  to  the  Boltzmann  distribution  before  the  reaction 
takes  place.  Furthermore,  lifetimes  of  the  excited  states  are 
smaller  the  higher  the  rovibrational  energy.  This  causes  the 
system  to  rapidly  return  to  the  ground  state.  The  expected 
effect  of  including  the  excited  states  in  the  dynamics  would 
be  a  slight  lowering  of  the  reaction  rate  coefficients.  There¬ 
fore,  in  this  work  we  first  consider  only  the  ground  state  PES 
since  it  dominates  the  reaction  dynamics  under  the  relevant 
physicochemical  conditions. 


The  smoothness  of  the  RKHS  interpolation  is  evident 
in  the  contour  plots  of  Figure  3.  In  panel  (a),  the  PES  for 
the  “NOl-l-02”  channel  is  shown.  The  global  minimum  was 
found  at  —  3.21  Oq  and  —  150.8°  (rj,jQ  =  2.25  ag  and 
^o-N-o  ~  133.1°  (Table  II)  in  internal  coordinates).  The  PES 
for  the  “0102-I-N”  channel  is  reported  in  panel  (b).  Two  sym¬ 
metric  minima  are  observed  for  the  same  Jacobi  distance 
which  is  measured  from  the  center  of  mass  of  O2  to  the  N 
atom. 

In  order  to  assess  the  quality  of  the  interpolated  RKHS 
PES,  an  additional  (off-grid)  cut  for  r  =  1.21  A  (N-Ol  dis¬ 
tance)  and  a  —  34°  was  calculated  along  the  R-coordinate  at 
the  MRCI-|-Q/cc-pVQZ  level.  The  comparison  between  the 
computed  ab  initio  energies  and  those  from  the  3d-interpolant 
is  reported  in  Figure  4.  The  average  mean  difference  over  the 
entire  range  is  0.03  eV  (see  inset  in  Figure  4).  The  cuts  for  two 
Gauss-Legendre  quadrature  angles,  that  is  (a  —  27.5  (open 
squares)  and  43.1°  (open  triangles)),  are  also  included  in  or¬ 
der  to  show  that  by  construction,  the  RKHS  curves  go  through 
the  ab  initio  energies. 
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TABLE  II.  Equilibrium  geometries  and  normal  mode  frequencies  from  the 
present  work  and  from  Ref.  13  and  experiments.  The  variables  rj  ^ 

6o_n-0’  equilibrium  N-0  distance,  the  O-N-0  angle,  and  the 

well  depth  for  the  global  minimum  of  NOj,  respectively. 


Equilibrium  NO2 

Present  work 

Ref.  13 

Experiments 

II 

>3’ 

0 

2.25 

2.26 

2.255*5 

^o-N-o  (deg) 

133.1 

134.3 

133.9*5 

D,  (eV) 

3.21 

3.11 

3  2314.51,52 

Normal  frequencies  (cm~^ 

) 

Symmetric  stretch 

1322.7 

1333.6 

1319.79450 

Asymmetric  stretch 

1623.3 

1658.0 

1616.85250 

Bending  mode 

759.6 

755.0 

749.64950 

Some  features 

of  the  present 

ground 

State  PES 

are  summarized  in 

Table  II  and 

compared 

to  results 

from  previous  calculations*^  (using  icMRCI+Q/cc-pVQZ 
with  CASSCF(13,10))  and  to  experiments*"*’ (high- 
resolution  microwave  spectroscopy*^  and  laser  induced 
fluorescence^*).  Good  agreement  with  previous  computations 
and  experiments  was  found  from  the  methods  used  in  the 
present  work. 

B.  Thermal  rate  coefficients 

A  total  of  10  000  individual  trajectories  was  run  at  each 
temperature,  including  T  —  100,  200,  300,  400,  1000,  1500, 
2000,  2700,  5000,  7000,  10  000,  15  000,  and  20  000  K.  In 
the  following,  the  thermal  rate  coefficients  obtained  from  the 
trajectory  calculations  are  presented.  The  convergence  of  the 
computed  rate  coefficient  is  reported  in  Figures  5(a)  and  5(b). 
In  panel  (a),  the  thermal  rate  coefficient  for  NO2  formation 
in  the  high-pressure  limit,  k°°(T),  is  shown  as  a  function  of 
the  number  of  trajectories  for  T  =  10000  K.  The  relative  er¬ 
ror  ((2cts/A:“)100%,  where  is  the  standard  deviation)  de- 


FIG.  4.  MRCI+Q/cc-pVQZ  energies  along  the  R  coordinate  for  r=  1.21  A 
(N-01  distance)  and  a  =  27.5°,  34.0°,  and  43.1°  (open  squares,  open  circles, 
and  open  triangles,  respectively).  R,  r,  and  a  are  Jacobi  coordinates.  The  solid 
lines  are  RKHS  interpolants.  The  inset  in  the  graph  represents  a  close-up  of 
the  cut  for  a  =  34.0°  which  is  not  a  Gauss-Legendre  quadrature  angle. 


FIG.  5.  Convergence  of  the  Monte  Carlo  integral  for  the  thermal  rate  coef¬ 
ficients  for  NO2  formation  in  the  high-pressure  limit  at  T  =  10000  K  versus 
the  number  of  trajectories  (panel  (a)).  The  confidence  interval  corresponds 
to  two  times  the  standard  deviation  (cr^).  In  panel  (b),  the  relative  error  as 
a  function  of  the  number  of  trajectories  is  reported  for  T  =  300  K  (solid), 
1000  K  (dashed  line),  5000  K  (dotted  line),  and  10  000  K  (dashed-dotted 
line).  The  red  horizontal  line  indicates  a  relative  error  of  10%. 

creases  to  ^2.5%  if  all  10  000  trajectories  are  used  to  compute 
k°°iT).  In  panel  (b),  the  decrease  of  the  relative  error  with  the 
number  of  trajectories  is  shown  for  selected  temperatures.  The 
horizontal  line  indicates  10%  of  the  relative  error  for  which 
R^hOO  trajectories  are  required.  This  indicates  that  the  com¬ 
putational  strategy  is  very  efficient  in  exploring  the  relevant 
regions  of  the  initial  phase- space  for  the  thermal  rate  coeffi¬ 
cient  calculations. 

In  Figure  6,  the  thermal  rate  coefficients  for  NO2 
formation  in  the  high-pressure  limit,  k°°{T),  are  shown  and 
compared  with  published  computational  and  experimental 
results.  Even  though  the  major  interest  here  concerns  the 
high-temperature  range,  these  data  provide  a  validation  of  the 
computational  procedure  put  forward  in  the  present  work.  We 
report  the  error  bar  as  twice  the  standard  deviation  (ct^)  in 
order  to  include  95.4%  of  the  cases  in  the  confidence  interval. 
In  the  following,  it  is  important  to  emphasize  that  the  present 
results  are  from  simulations  using  the  l^A'  PES  only  (open 
squares  with  solid  line).  The  present  results  are  in  reasonable 
agreement  with  results  from  previous  calculations*^  (solid 
down-oriented  triangles)  which  also  used  the  l^A'  PES  only, 
although  with  a  different  parametrization.  Both  calculations 
underestimate  the  experimentally  determined  thermal  rate 
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FIG.  6.  Thermal  rate  coefficients  for  the  NOj  formation  in  the  high-pressure 
limit,  k^{T)  versus  temperature.  In  the  legend,  [a]  =  Ref.  17,  [b]  =  Ref.  16, 
[c]  =  Ref.  18,  [d]  =  Ref  53,  [e]  =  Ref  54,  and  [f]  =  Ref  58. 

coefficients  at  T  —  200,  300,  and  400  K  (solid  circles  and 
squares).  At  2700  K,  the  computed  rate  coefficient  compares 
favourably  with  experimental  data.^^’^"^ 

No  correction  was  made  for  zero-point  energy  (ZPE)  and 
other  quantum  effects,  such  as  tunneling  through  the  centrifu¬ 
gal  barrier  which  can  be  still  significant  at  low  temperatures. 
The  centrifugal  barrier  originates  from  the  conservation  of  the 
total  angular  momentum  and  consists  of  the  rotational  energy 
associated  with  the  orbital  angular  momentum  of  the  colliding 
partner  with  respect  to  the  center  of  mass  of  the  whole  sys¬ 
tem.  The  orbital  angular  momentum  is  directly  proportional 
to  the  impact  parameter  and  with  the  initial  relative  momen¬ 
tum  of  reactants.  An  increase  of  the  thermal  rate  coefficient 
for  lower  temperatures  has  been  found  to  be  related  to  tun¬ 
neling  through  the  entrance  barrier  in  the  reaction  of  OH  with 
methanol  (T <  200  K)  and  0(^P)  with  alkenes  (T-  20 
This  has  also  been  shown  in  theoretical  studies  of  HCOH 
isomerizations.^^  On  the  other  hand,  results  when  using  the 
two  lowest  l^A'  and  l^A"  PESs  indicate  that  the  l^A"  PES 
might  play  an  important  role  along  with  quantum  effects  at 
temperatures  below  400  K.  However,  for  temperatures  above 
1000  K  the  results  using  only  the  l^A'  PES  reasonably  match 
with  the  experimental  results,  e.g.,  for  1500  and  2700  K.  This 
can  be  explained  by  noting  that  the  well  of  the  electronically 
excited  C-state  of  NO2  (C  l^A")  is  comparatively  small  and 
leads  to  redissociation  of  NO2  at  high  temperatures.'^  Given 
this,  the  influence  of  tunneling  and  the  l^A"  PES  at  high  tem¬ 
perature  appear  to  be  negligible. 

For  practical  work  in  reaction  networks,*^  it  is  advanta¬ 
geous  to  use  a  parametrized  form  for  ld^{T).  A  previously  con¬ 
sidered  empirical  expression'^  has  been  fitted  to  the  results 
from  simulations  using  the  l^A'  PES, 

yfc“(r)/(10“"cmVs) 

=  2.715 -f  0.4490  In  T/K-f  0.3385  (In  T/Kf 
-0.7134  (In  r/K)l  (17) 

This  yields  the  red  solid  line  in  Figure  6. 


FIG.  7.  High-pressure  limit  thermal  rate  coefficients  for  NO2  formation 
(squares)  (fit  with  black  solid  line),  NO2  formation  and  dissociation  (cir¬ 
cles),  oxygen  exchange  channel  (up-oriented  triangles),  and  oxygen  produc¬ 
tion  channel  (diamonds). 


FIG.  8.  Thermal  reaction  rate  coefficients  for  NO2  formation  and  decay, 
for  different  thresholds  of  the  lifetime.  Results  for  lower  tempera¬ 
tures  (100-2700  K)  are  shown  in  panel  (a)  and  for  higher  temperatures  (5000- 
20  000  K)  in  panel  (b). 
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FIG.  9.  (a)  and  (b)  Ratios  of  the  rate  coefficients  of  different  channels. 


The  thermal  rate  coefficients  for  temperatures  between 
T  —  5000  K  and  T  —  20  000  K  are  reported  in  Figure  7  for 
cases  TIV  studied  in  this  work.  With  increasing  temperature 
a  slight  increase  in  the  thermal  rate  coefficients  is  found.  At 
20  000  K,  they  increase  to  =  5.5  x  10“"  cm^/s, 

=  3.4  X  10“^'  cm^/s,  =  1.9  x  10”'*  cm^/s,  and  ^^2+^^ 
=  0.65  X  10”"  cm^/s.  Again,  for  practical  purposes  the  data 
for  k°°iT)  were  fitted  to  a  linear  equation  as  (solid  line  in 
Figure  7) 

k“(r)/(10”"cmVs)  =  2.497 +  0.1529(r/K).  (18) 

Figure  8  reports  the  rate  coefficient  for  NO2  formation 
and  decay,  k^^2{T),  for  different  thresholds  for  the  lifetime 
of  the  NO2  molecule.  Panel  (a)  shows  results  for  lower  tem¬ 
peratures  (100-2700  K)  and  panel  (b)  those  for  higher  tem¬ 
peratures  (5000-20  000  K).  The  results  indicate  that  the  NO2 
molecule  lifetimes  (see  definition  of  lifetime  (f^fe)  in  Sec.  II) 
are  mostly  below  16  ps  and  1.6  ps  for  low  and  high  tempera¬ 
tures,  respectively. 

In  order  to  identify  which  process  is  relevant  at  higher 
temperatures,  the  ratios  between  the  thermal  rate  coefficients 
are  shown  in  Figures  9(a)  and  9(b)  for  the  different  cases. 
In  all  cases,  the  ratios  become  almost  constant  at  higher 
temperatures.  The  values  for  20  000  K  are:  k°°lkf^^  Rs  3, 
j^NO,/^ex  ^  ^  1.6,  k“/k°2+N  sa  8.4, 

k^^2/k^2+^  Rs  5.2,  and  k“/k°2+N  sa  2.8.  At  this  temperature. 


FIG.  10.  Fraction  of  initial  vibrational  energy  of  NO  with  respect  to  the  total 
energy  versus  the  vibrational  quantum  number  for  8  values  of  the  collision 
energy  (panel  (a)).  The  contour  plot  of  the  channel  probability  as  percentage 
is  shown  in  panel  (b).  Only  data  for  “O2+N”  channel  are  shown. 


the  oxygen  exchange  channel  is  about  3  times  faster  than  the 
“O2-I-N”  channel. 

For  a  more  detailed  description  of  the  O2  production 
channel.  Figure  10  (panel  (a))  reports  the  fraction  of  initial 
vibrational  energy  of  NO  relative  to  the  total  energy  versus 
the  vibrational  quantum  number  of  NO  (%o)  ^  values  of 

the  collision  energy.  The  fraction  of  collision  energy  is  com¬ 
plementary  to  the  vibrational  energy  to  give  the  total  energy 
of  the  system  +  ^vib  —  ^tot)'  In  panel  (b),  the  contour 
plot  of  the  channel  probability  is  shown.  In  all  cases,  the  rota¬ 
tional  ground  state  —  0)  is  assumed.  For  this  calculation, 
2500  trajectories  were  run  for  each  combination  of  ={0.82, 
1.09,  1.36,  2.18,  2.72,  4.08,  4.63,  5.44}  eV  and  ={0’  1’ 
2,  3,  4,  6,  9,  12,  15,  20}  values.  The  results  indicate  that  for 

=  0  the  02-channel  opens  only  for  very  high  collision 
energies  (>4.63  eV).  However,  for  an  initial  NO-state  with 

=  9  (corresponding  to  E^  =  0.82  eV)  production  of  0.6% 
of  O2  is  observed.  This  qualitatively  follows  from  Polanyi’s 
rule'"  which  predicts  that  the  reaction  should  go  to  the  late- 
barrier  “O2  -t-  N”  channel  when  the  NO  molecule  is  initially 
vibrationally  highly  excited  even  for  low  collisional  energies. 
In  order  to  give  a  clearer  picture  of  this,  the  reaction  is  con¬ 
sidered  to  start  from  N-I-O2  (left  side  in  Figure  1).  In  this 
situation,  the  reactants  experience  an  early  barrier  of  about 
1.3  eV  (Figures  1  and  3(b)).  The  exit  channel  in  this  case 
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is  barrierless  (NOl  +  02).  According  to  Polanyi’s  rule,  the 
excess  energy  released  after  crossing  the  early  barrier  goes 
predominantly  into  the  internal  degrees  of  freedom  in  the 
products  (NOl)  rather  than  into  their  (relative)  translational 
energy. 

For  the  set  of  vibrational  quantum  numbers  and  ener¬ 
gies  studied,  the  maximum  channel  probability  of  13%  for 
02-formation  is  reached  at  =  20.  In  general,  the  trend 
indicates  that  for  lower  vibrational  quantum  numbers,  more 
energy  needs  to  be  in  the  translational  degree  of  freedom. 
Conversely,  for  higher  vibrational  excitation  of  the  NO,  less 
collision  energy  is  required.  When  both  the  collision  energy 
and  the  vibrational  quantum  number  of  NO  are  high,  the  re¬ 
action  probability  of  N-f02  formation  is  maximal. 

IV.  CONCLUSIONS 

The  0(^P)  -t-  NO(^n)  reaction  has  been  studied  over 
a  range  of  temperatures  relevant  for  the  hypersonic  flight 
regime  (5000-20000  K).  For  this  purpose,  the  thermal  cross 
sections  and  rate  coefficients  were  calculated  for  the  dif¬ 
ferent  exit  channels.  The  results  indicate  that  oxygen  ex¬ 
change  (N02-I-01)  occurs  at  all  temperatures.  However,  02- 
formation  is  only  found  for  T  >  5000  K.  For  both  channels, 
the  thermal  rate  coefficient  increases  as  the  temperature  in¬ 
creases  {T  —  5000-20  000  K).  A  more  detailed  scrutiny  of 
the  “0102-l-N”  channel  indicates  that  the  smaller  the  vibra¬ 
tional  quantum  number  of  NO,  the  higher  translational  energy 
is  required  to  open  it;  conversely,  the  higher  the  vibrational 
quantum  number,  the  less  translational  energy  is  needed.  Our 
results  for  lower  temperatures  are  discussed  and  compared  to 
previous  experimental  and  computational  results.  In  this  re¬ 
spect,  this  work  also  presents  the  validation  of  an  efficient 
computational  protocol  for  the  investigation  of  reaction  dy¬ 
namics  of  triatomic  systems.  Two  aspects  make  the  approach 
efficient  and  accurate.  First,  the  ab  initio  energies  are  handled 
by  means  of  a  RKHS  representation  which  is  accurate,  com¬ 
putationally  convenient,  and  correctly  describes  the  asymp¬ 
totic  regions  of  the  PESs.  Second,  thermal  quantities  are  cal¬ 
culated  by  using  Importance  Sampling  which  is  efficient  to 
evaluate  multi-dimensional  integrals  over  the  phase-space  of 
initial  states.  With  less  than  1000  classical  trajectories  a  rela¬ 
tive  error  below  10%  for  the  rate  coefficients  is  obtained. 
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